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ABSTRACT 


This  paper  presents  the  numerical  calculation  of  the  Wannier  functions  (an(z))  of 
the  GaAsZAlo.25Gao.75As  superlattice  structure.  The  Wannier  functions  are  linear 
combinations  of  the  Bloch  functions  ( y/^  (z) )  that  can  be  viewed  as  a  convenient 

mathematical  instrument  to  get  around  the  lack  of  orthogonality  of  the  tight  binding 
formulation  and  are  useful  tool  when  the  position  of  an  electron  has  physical  importance. 
However,  except  for  finding  a  variational  principal  for  the  energy  levels  of  one 
dimensional  crystals  in  terms  of  The  Wannier  functions  and  simple  cubic  lattice 
calculations  for  instructional  reasons,  models  have  not  been  calculated  or  plotted  showing 
the  theory  can  be  expanded  to  more  complex  superlattice  structures.  We  develop  an 
algorithm  to  calculate  numerically  and  plot  the  Wannier  functions  of  the 
GaAsZAlo.25Gao.75As  superlattice  structure  and  hence  prove  that  these  functions  can  be 
calculated  even  for  complex  structures.  By  using  the  plots  of  our  numerical  modelling, 
we  displayed  the  peculiar  properties  of  the  Wannier  functions: 

a.  Wannier  functions  are  real. 

b.  They  fall  off  exponentially. 

c.  They  are  either  symmetric  or  anti-symmetric  about  z- 0. 


VI 


TABLE  OF  CONTENTS 

l.  INTRODUCTION . 1 

n.  THE  W ANNIER  FUNCTIONS  AND  THEIR  PROPERTIES . 3 

m.  SUPERLATTICE  EFFECTIVE-MASS  SCHRODINGER  EQUATION . 7 

IV.  THE  WANNIER  FUNCTIONS  IN  FINITE  DIFFERENCE  FORM . 13 

V.  RESULTS  AND  CONCLUSIONS . 17 

A.  BANDSTRUCTURE  FUNCTION . 17 

B.  WANNIER  FUNCTIONS . 19 

C.  CONCLUSION  AND  FINAL  REMARKS .  24 

APPENDIX  A:  COMPUTER  CODE  IN  MATLAB  PROGRAMMING  LANGUAGE 
FOR  NUMERICAL  CALCULATION  OF  THE  BANDSTRUCTURE 
FUNCTION . 29 

APPENDIX  B:  COMPUTER  CODE  IN  MATLAB  PROGRAMMING  LANGUAGE 
FOR  NUMERICAL  CALCULATION  OF  THE  WANNIER 
FUNCTIONS . 33 

LIST  OF  REFERENCES . . .  37 

INITIAL  DISTRIBUTION  LIST . . 39 


vii 


viii 


ACKNOWLEDGEMENT 


First  of  all,  I  would  like  to  thank  my  advisor  Prof.  J.  H.  Luscombe  who  made  this 
lifetime  goal  come  true.  Secondly,  I  would  like  to  express  my  gratitude  to  Prof.  R.  L. 
Armstead  for  his  support. 

I  would  like  to  dedicate  this  thesis  to  my  family,  my  wife  Birgul  and  my  son 
Egem,  without  whom  this  goal  cannot  be  achieved.  She  has  showed  me  that  the  biggest 
input  to  success  is  a  loving  family. 

And  finally  I  want  to  thank  my  good  friend  Scott  Bewley  (USN)  and  his  wife 
Lorelle  who  showed  me  what  a  good  friendship  is  and  who  helped  me  and  supported  me 
in  all  phases  of  this  thesis. 


IX 


X 


I.  INTRODUCTION 


Semiconductor  superlattices  are  a  novel  class  of  man-made  solids  which  consist 
of  a  repetition  of  alternating  layers  of  two  or  more  different  semiconductor  materials.  By 
adding  artificial  periodicity  to  the  multi-layered  structure  [Ref.  1],  the  electronic  nature  of 
the  superlattice  can  be  controlled  through  the  constituent  material’s  crystalline 
orientation,  layer  thickness  and  material  composition.  The  result  of  this  periodic  layering 
process  creates  semiconductor  superlattices  which  present  exciting  optical  and  electronic 
characteristics  which  cannot  be  found  in  existing  semiconductor  crystals  [Ref.  2]. 

From  the  beginning  of  the  study  of  this  new  class  of  solids,  most  of  the  attention 
and  research  has  been  concentrated  on  new  generations  of  optoelectronic  device 
applications.  One  of  the  promising  device  implementations  is  fast  and  highly  sensitive 
InAlGaAs/InAlAs  superlattice  avalanche  photodiodes  for  use  in  1.3  to  1.6  pm  optical 
communications  [Ref.  3].  Other  areas  of  research  include  the  study  of  superlattices  for 
the  future  generation  of  far  infrared  detector  imaging  arrays  with  wavelength  X>12  pm 
[Ref.  4]  and  the  study  of  high  optical  confinement  effects  [Ref.  5]  leading  to 
investigations  of  using  superlattices  to  reduce  the  threshold  current  density  and  enhance 
the  device  stability  of  semiconductor  lasers. 

Even  though  successful  device  implementations  have  been  reported,  and  a  world¬ 
wide  research  effort  is  underway  in  this  new  field  of  study,  many  problems  remain 
unsolved.  These  problems  include  everything  from  the  preparation  of  the  superlattices  to 
understanding  the  transport  and  optical  phenomena  introduced  by  these  novel  materials. 

In  this  thesis,  the  Wannier  functions  (a„  ( z-dj ))  of  a  GaAsZGao.75Alo.25As 
superlattice  is  numerically  calculated  and  plotted  with  the  new  lattice  constant  in  the 
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growth  direction  z  set  by  ^GaAs+^AiGaAs  (Fig-  !)•  The  Wannier  functions  can  be  viewed  as 


V(z),  eV 


Figure  1 :  Conduction  band  edge  profile  of  GaAs/Alo.75Gao.25As  superlattice  with  the  new 
lattice  constant  d  =  ^AiGaAs  +  ^GaAs  where  c^AiGaAs  an(i  ^GaAs  are  are  thickness  of  the 
corresoondine  semiconductors. 

a  convenient  mathematical  instmment  to  get  around  the  lack  of  orthogonality  of  the  tight- 
binding  formulation .  It  is  also  a  useful  tool  in  understanding  the  spatially  localized 
nature  of  electrons  when  they  are  near  an  impurity  state.  Except  for  theoretical  studies  on 
finding  a  variational  principal  for  the  energy  levels  of  one  dimensional  crystals  in  terms 
of  the  Wannier  functions  [Ref.  6],  and  simple  cubic  lattice  calculations  for  instructional 
reasons  [Ref.  7],  models  have  not  yet  been  calculated  or  plotted  showing  the  theory  can 
be  expanded  to  more  complex  superlattice  designs.  This  thesis  will  numerically  calculate 
and  plot  the  Wannier  functions  of  GaAs/Al0.25  Gao.75  As  superlattice  structure  to  prove  the 
theory  beyond  the  simple  cubic  problem  and  will  display  the  properties  of  these 
functions  graphically. 
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n.  THE  WANNIER  FUNCTIONS  AND  THEIR  PROPERTIES 

In  a  one  dimensional  crystal,  the  mathematical  basis  for  the  motion  of  an  electron 
can  be  expressed  by  the  Schrodinger  equation  as 

ft 2  d2 

- —j+V(z  +  c)  \jr  =  Ey  (1) 

2m  dz 

where  V  (z)  is  a  periodic  potential  function  of  the  crystal  in  the  growth  direction  z  with 
the  period  of  lattice  constant  c,  i.e.,  V(z)=V(z+c)  and  m  is  the  mass  of  the  electron.  Bloch 
has  shown  that  the  eigenfunctions  of  the  one-electron  Schrodinger  equation,  or  Bloch 
states,  for  this  periodic  potential  (Eqn.  1)  have  the  mathematical  form  [Ref.  8] 

¥,±(z)  =  eikzunk(z)  (2) 

where  «nk(z)  is  a  function  with  the  periodicity  of  the  potential.  The  subscript  n  stands  for 
the  n*  energy  band  to  which  the  eigenfunctions  belong  and  the  subscript  k  is  the  wave 
vector  in  the  growth  direction  z. 

In  the  method  of  the  tight  binding  formulation,  using  atomic  wave  functions  as  a 
basis  set  has  a  practical  advantage  over  Bloch  states  in  qualitative  calculations  [Ref.  9]. 
The  actual  calculations  with  those  functions,  however,  create  difficulties  due  to  the 
atomic  wave  functions  being  centered  on  different  atomic  sites  that  are  not  orthogonal  to 
each  other.  Furthermore,  they  are  not  necessarily  orthogonal  for  different  k  values  either. 

Wannier  has  introduced  a  set  of  localized  wave  functions,  an{z-cj),  which  are 
linear  combinations  of  all  the  Bloch  states  (z)  for  a  given  energy  band,  and  have  the 
advantages  of  orthogonality  [Ref.  10]: 

an(z-Cj)  =  -J=  £  e~lkCi  y/^  (z)  (3) 
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where 


a"^  =  VF  ^ 

N  is  the  number  of  lattice  sites  in  the  one  dimensional  crystal  and  subscript;  stands  for 
the  lattice  site  of  which  we  are  going  to  calculate  the  Wannier  function.  The  summation 
in  Equation  4  extends  over  all  values  of  the  first  Brillouin  zone  (BZ).  The  Bloch  function 
\f/nk  (4)  is  continuous  in  k,  and  is  the  same  function  of  k  at  k+K  where  K  is  a  reciprocal 

lattice  vector,  i.e.,  the  Bloch  functions  are  periodic  in  the  reciprocal  lattice  [Ref.  11]. 
Since  the  Wannier  functions  are  functions  of  distance  from  the  direct  lattice  site;,  we 
may  think  of  Equation  3  and  4,  i.e.,  y/nk  (z)and  an(z-cj)  as  the  Fourier  transform  of  each 
other.  The  Fourier-transform  of  Equation  3  defining  the  Wannier  functions,  can  now  be 
calculated  to  obtain  the  Bloch  functions: 

<5> 

The  Wannier  functions  can  show  orthonormality  properties  in  contrast  to  the  tight 
binding  formulation,  because  the  Wannier  functions  centered  on  different  lattice  site 
and/or  or  belonging  to  different  energy  bands  are  orthogonal  to  each  other : 

J a* (z - Cj )am (z - c , )dz  =  -j=^e,kCi e~,kCi J yrA  (zty^- ( z)dz 

=  J-Y  eiikCrkCl)S  S  ■ 

Nh 


(6) 


where  we  have  invoked  the  summation  relation 

JV"'-*'1  =  AW  (7) 

k 

and  the  orthogonality  of  the  Bloch  functions.  Hence,  the  Wannier  functions  form  a 
complete  orthonormal  set,  and  the  completeness  can  be  shown  by  summing  over  all 

energy  bands  and  all  lattice  sites, 

£  al  (z  ~  cj K  (z ’  ~  cj )  =  S{z  -  z ') .  (8) 

For  a  one  dimensional  periodic  potential  with  no  degeneracies  of  energy  bands,  It 
has  been  rigorously  shown  that  the  Wannier  functions  will  decay  exponentially  at  large 
distances  [Ref.  12].  For  our  purposes,  we  can  easily  show  this  localized  nature  of  the 
Wannier  functions  qualitatively.  If  we  replace  z  with  c,  in  Equation  3  where  c,  is  a  lattice 
vector,  then  we  get 

an (c,  - C; )  =^Yie~'KcrCl)ynk  (z)  (9) 

When  I  Cj  -  Cj  I  becomes  large  the  exponential  part  of  the  summation  oscillates  rapidly  and 
contributes  vanishingly  to  the  summation.  Hence,  we  can  expect  the  Wannier  functions 
of  one  energy  band  to  fall  off  exponentially  at  large  distances  and  to  be  identical 
functions  which  are  centered  at  each  of  the  lattice  points. 

One  of  the  most  interesting  properties  of  the  Wannier  functions  is  their  reality. 
Regardless  of  the  phases  of  the  Bloch  functions,  Wannier  functions  are  real.  W.  Kohn 
has  discussed  the  reality  problem  of  the  W annier  functions  in  detail  [Ref.  12].  We  are  not 
going  to  prove  this  peculiar  aspect  of  the  Wannier  functions  in  this  paper,  but  we  refer  the 
reader  to  W.  Kohn’s  paper. 
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With  the  aforementioned  properties,  the  Wannier  functions  can  be  considered  a 
useful  mathematical  tool  when  the  position  of  an  electron  has  physical  importance,  e.g., 
in  the  discussion  of  the  motion  of  an  electron  around  an  impurity  site.  Other  important 
areas  of  application  include  deriving  a  transport  theory  for  Bloch  electrons  [Ref.  13]. 
Also,  using  Wannier  functions  in  the  calculation  of  energy  bands  has  been  discussed  but 
no  practical  application  emerged  [Ref.  14]. 
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in.  SUPERLATTICE  EFFECTIVE-MASS  SCHRODINGER 

EQUATION 

Before  presenting  our  method  of  numerical  calculation  of  the  Wannier  functions, 
we  are  going  to  derive  the  Schrodinger  equation  for  the  GaAs/Alo.7sGao.25As  superlatice 
structure  (Fig.  1)  by  using  the  effective  mass  theory  with  the  Kronig-Penney  model.  This 
approach  describes  the  electron  with  an  effective  mass  m*  where  subscript  i  corresponds 

to  the  constituent  layers  of  the  superlattice,  and  approximates  the  superimposed 
perturbing  potential  (Superlattice  bandedge  profile)  as  periodic  square  wells  (Fig.  1). 

If  H  is  the  Hamiltonian  for  an  electron  in  a  one  dimensional  perfect  crystal,  then 
H¥nk(z)=En(k)¥nk(z )  (10) 

where 

H  =  ~^y+V(z),  (11) 

2m  dz 

and  En(k)  are  the  eigenvalues  of  the  Bloch  functions  ^(z),  namely  the  band  structure 
function.  Since  band  structure  function  EJk)  is  a  non-degenerate  and  continuous 

function  of  k  with  the  periodicity  of  the  crystal  in  the  reciprocal  space  (fc-space)  [Ref.  15], 
we  have  the  possibility  of  defining  the  Fourier  expansion  of  it  as  we  did  for  the  Bloch 
functions  (Eqn.  5): 

£„(*)  =  £  E„(c)e“  •  '  (12) 

where  the  summation  is  over  the  lattice  sites  of  the  host  material  in  the  growth  direction 
(z),  and  En  (c)  is  the  Fourier  expansion  coefficient. 
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Since  the  band  structure  function  En(k)  is  an  analytic  function  of  k  and  the  wave 
functions  y/nk  (z)  have  the  periodicity  of  the  crystal,  we  can  assume  the  replacement  of 

0 

the  perfect  crystal  Hamiltonion  (H)  by  the  operator  En  {—i  — )  which  follows  from 


o 

En(k )  when  Hermitian  differential  operator  —  i—  is  used  in  place  of  k  [Ref.  16]  where  i 


is  V-l  •  We  will  show  that  H  and  E  (-1— ) 

oz 


have  the  same  eigenvalues  and 


E(-i ^-)Wnk  (z)  =  'ZEn(C)'/'nk(Z  +  C) 

wZ  c 

c 

=  En(k)Wnk(z)  .  (13) 

We  have  used  the  Fourier  transform  of  the  band  structure  function  (Eqn.  12)  and  Bloch’s 
theorem  [Ref.  17]: 

¥ttk(z±c)  =  e±*cynk(z)  .  (14) 
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As  can  be  noted,  the  periodic  potential  of  the  crystal  structure  V(z)  in  Equation  1 1  is  not 
explicitly  present  in  this  formulation  (Eqn.  13)  where  it  is  incorporated  with  the  usual 
kinetic  energy  operator, 


h2  a2 


(15) 


2m  dz2 

The  unperturbed  Hamiltonian  defined  by  the  modified  kinetic  energy  operator 
'  h2  d2 


E(-i — )  would  be  far  more  complicated  than  the  kinetic  energy  operator  -- 


,  if 


dz'  •”  1  w'  *  2  mdz2 

we  tried  to  include  the  full  ^-dependence  of  the  band  structure  function  En(k)  [Ref.  18]. 

However,  many  aspects  of  the  electronic  properties  of  semiconductors  can  be  obtained 
from  a  small  number  of  excited  electrons  found  at  the  bottom  of  the  conduction  band  and 
holes  created  in  the  absence  of  these  electrons  at  the  top  of  the  valence  band.  Hence,  what 
we  need  is  the  form  of  the  band  structure  function  for  these  carriers  over  a  small  region  of 
fc-space  around  the  band  extrema.  This  implies  that  the  bandstructure  function  can 
generally  be  drawn  as  a  quadratic  function  of  k.  Besides,  for  most  of  the  semiconductors, 
En(k)  is  approximately  parabolic  near  an  extremum,  as  for  the  conduction  band  edge  of 
GaAs  and  Alo.75Gao.25 As  [Ref.  19].  For  the  present  problem,  by  restricting  our  treatment 
to  the  band  extremum  E0  occuring  around  the  wave  vector  value  k  =  k0  =  0  which  is  the 
case  for  direct  gap  semiconductors  like  GaAs,  we  take  the  Taylor  series  expansion  of 

EM- 


E„(k)  =  E0  +  P,(k-  k0)  +  ^fi2(k-  k0)2  + (16) 

where  the  expansion  coefficients  are 
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(17) 


Pi  = 


(d‘E 

dk‘ 


Jk~ko 


the  subscript  i  representing  the  order  of  expansion.  Since  the  Taylor  expansion  of  E-k 
relation  around  the  band  extremum  is  relevant  to  the  problem,  the  first  order  term 
(k-k0)  will  vanish  by  construction,  and  higher  order  terns  can  be  neglected.  For 

dimensional  reasons,  (32  must  have  the  units  of  ft2  / m .  Therefore,  by  defining  an 
effective  mass  (m*)  of  the  electron  and  replacing  the  curvature  of  the  E-k  relation  with 


ft2  /  2 m  ,  the  band  structure  function  near  the  extremum  point  can  be  characterized  as 

£„(*)  =  E0+^*!  .  (18) 

2m 


Once  again,  by  replacing  the  canonical  substitution  -  i—  in  place  of  k  in  Equation  18, 
we  obtain  the  effective  mass  Schrodinger  equation  for  a  bulk  crystal. 


-a2  a2 

2 m  dz2 


+  E0 


¥nk(z)  =  En(k)y/nk(z) 


(19) 


where  the  effect  of  the  crystal  potential  V(z)  is  absorbed  in  the  parameter  m*. 

Now  that  the  Schrodinger  equation  for  the  perfect  crystal  has  been  formulated,  we 
can  introduce  our  perturbing  potential  Us(z)  shaped  by  the  GaAs/Alo.25Gao.75As 
superlattice  in  the  growth  direction  z.  By  applying  the  Kronig-Penney  model  depicted 
earlier,  the  superimposed  superlattice  potential  can  be  expressed  as 
U,  (z)  =  EfGaAs  -  EGaAs  =  A Ec  \z-  nd\)dGaAs  /  2 

=  0  \z-nd\(dGaAs/2  ,  (20) 
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where  EGaAs  and  E*lGaAs  are  the  bottom  of  the  conduction  bands  of  GaAs  and 
Alo.25Gao.75As  respectively  and  d  is  the  superlattice  constant  defined  earlier.  Equation  20 
explicitly  implies  that  Us(z )  is  aperiodic  function  of  z  with  the  periodicity  d,  i.e., 

Us(z)  =  Us(z  +  nd)  .  (21) 

By  defining  the  perturbing  superlattice  periodic  potential  Us  (z) ,  which  is 
superposed  over  the  intrinsic  crystal  potential  of  the  host  material,  the  Schrodinger 
equation  for  GaAs/Alo.25Gao.75As  superlattice  can  be  then  written  as 

- Stt It ti'.ft)  (*><",» <z) ,  (22) 

2m  (z)  dz 

where  we  moved  E0  to  the  right-hand  side  of  the  Equation  22  and  defined  £n  (k)  which 
is  the  difference  between  bandstructure  function  En  (k)  and  band  extremum  Eg .  This 


rearranging  of  energy  terms  will  give  the  energy  values  from  the  edge  of  the  conduction 
band  of  the  constituent  layer  (Since  EGaAs  <  EfGaAs ,  in  this  context  it  is  GaAs)  which  is 
approximated  as  a  quantum  well  by  the  Kronig-Penney  model.  The  spatial  dependence  on 
the  growth  direction  z  is  brought  in  to  the  effective  mass  (m*(z))  since  the  materials 
forming  the  layers  of  the  superlattice  have  different  values  of  effective  mass. 

Since  the  Bloch  functions  are  orthogonal  for  different  energy  values  which  are 
real,  it  requires  the  hermiticity  of  the  Hamiltonian  (Energy  operator). 


*2  d2 

H  = - " — ^T+C/i(z) 

2m* (z)  dz2 


But,  this  form  of  the  Hamiltonian  is  not  accurate  when  the  effective  mass  has  spatial 

dependence  ,  since  the  canonical  operator  —  i—  does  not  commute  with  z,  which 

dz 
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9  1 

requires  the  correct  ordering  of  —  i—  and  — ; —  terms.  For  example,  let  us  consider  an 

dz  m  ( z ) 

interface  between  the  two  constituent  layers  (GaAs  and  Alo.25Gao.75As),  say  at  z0  ■  The 
wave  function  ynk{z)  has  to  obey  the  usual  boundary  conditions  [Ref.  20]  to  guarantee 
current  conservation  across  the  junction, 


vMaAs)=vAzraAs) 


AlGaAs  \ 


(24) 


dz 


dVnk(z) 


ZGaAs 
0 


dz 


(25) 


y  AlGaAs 
c0 


where  zGaAs  and  z£lGaAs  are  the  sides  of  the  junction  in  GaAs  and  Alo.25Gao.75  As  layers 

respectively.  On  the  other  hand,  the  layers  have  different  effective  masses  [Ref.  21].  So, 
the  boundary  condition  for  matching  the  derivatives  of  the  wave  functions  (Eqn.  25)  must 
take  the  effective  masses  into  account: 


1 


* 


m 


GaAs 


dz 


ZGaAs 
0 


1 

* 

m  AlGaAs 


dz 


Z  AlGaAs 
0 


(26) 


where  m*GaAs  and  m\lGaAs  are  effective  masses  of  the  electron  in  the  corresponding 

semiconductors.  Since  the  derivative  in  Equation  26  is  apparently  the  Hermitian 
differential  operator,  the  hermiticity  here  also  guarantees  the  current  conservation.  As  a 
result  of  these  consideration,  the  correct  form  of  the  Schrodinger  equation  with  a  spatially 
dependent  effective  mass  is  given  by 

-»*j-  ‘  =  .  (27) 

dz  m  (z)  dz 

which  assures  the  hermiticity  and  guarantees  the  conservation  of  current  at  the  interface. 
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IV.  THE  WANNIER  FUNCTIONS  IN  FINITE  DIFFERENCE  FORM 


In  this  chapter  we  will  present  the  numerical  approach  which  is  utilized  to  obtain 
the  Wannier  functions  of  the  GaAs/Al0.25Gao.75As  superlattice  structure.  The  algorithm  is 
then  built  in  two  steps:  The  first  step  is  to  solve  Equation  27  for  the  wave  functions 
\ffnk  ( z )  and  the  eigenvalues  £n(k)  (Bandstructure  function).  In  the  second  step  we  only 

need  to  use  Equation  4  to  calculate  and  plot  the  Wannier  functions;  i.e.,  for  each  band 
index  n,  we  sum  the  wave  functions  over  all  values  of  k. 

We  are  going  to  use  finite  difference  approximation  [Ref.  22]  to  solve  the 
superlattice  effective  mass  equation  (Eqn.  27)  for  the  wave  functions.  In  this  method,  we 
evaluate  the  values  of  y/nk  (z)  at  discrete  points  Zj  which  are  obtained  by  dividing  the 


superlattice  unit  cell  into  finite  number,  N-l,  of  intervals  Az  defined  by 

Az  =  zj+l-Zj  . 


(28) 


The  first  derivative  in  the  finite  difference  form  at  the  point  zj  may  be  written  as 

3  =r^+Az)-MzcA*>  (2 

az^*^  2&z 

The  expression  for  the  second  derivative  can  be  obtained  by  applying  Equation  29 
successively  [Ref.  23]: 


3  1  a 

dz  m’(z)  dz 


Zj 


1 


m  (zj+l) 


VnkiZi+l)- 


m*  (zj+ 1)  m\zH) 


K*(Z;)  + 


m\zh ,) 


(2Az) 


(30) 


Before  writing  the  eigenvalue  matrix  equation  for  Equation  27,  we  have  to  take 
care  of  the  boundary  conditions  along  the  growth  direction  z.  In  the  artificial  periodicity  d 
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of  the  superlattice,  the  wave  function  must  satisfy  the  Bloch  theorem  (Eqn.  14)  [Ref.  17]: 

Wnk(z±d)  =  e±ikd¥nk(z)  .  (14) 

In  the  finite  difference  matrix  equation,  we  replace  ysnk(z0)  and  ynk{zN+x)  according  to 
Equation  14.  Since  we  have  expressions  for  each  mesh  point  Zj  in  the  superlattice  unit 
cell,  Equation  27  can  be  rewritten  as  an  eigenvalue  equation.  By  setting 
2Az  =  CO  =  d  /( N  - 1)  and  an  replacing  the  constants  in  Equation  27  with 


r  _  -  r  r  =. 

1  AlGaAs  ~  ~  *  7  A  1  GaAs 

^MGaAs® 


'^mGaAs® 


T  =  ro 


VAlGaAs  =&EC+- 


V  =■ 

Y  GaAs 


and  using  Equation  30,  we  get  the  eigenvalue  matrix  equation  for  each  value  of  k: 


~2Ta+Va 


0  0 


¥nk(Z2) 


■2Ta+Va 


¥nk(Z l) 

Wnk(z) 


0  rG  - 
o  o 
o  o 
o  o 
o  o 


2rG+vG 


2rG  +vG 


0  0 

0  0 

0  0 

0  0 

0  .0 


o  r\ 


~2Ta+Va 


¥nk(ZN) 


¥nk(ZN). 
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Once  we  have  obtained  the  wave  functions  for  each  value  of  the  reciprocal  wave 
vector  k  for  a  specific  band  index  «,  we  will  create  a  subroutine  to  sum  the  wave 
functions  over  all  the  k  values  in  the  first  Brillouine  zone  (BZ)  (Eqn.  4)  to  calculate  and 
plot  the  Wannier  functions  of  our  superlattice  structure. 
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V.  RESULTS  AND  CONCLUSIONS 

The  computer  codes  in  the  Appendix  A  and  B  have  been  developed  for  the 
solution  of  the  effective  mass  Schrodinger  equation  for  the  GaAs/Alo.25Gao.75As 
superlattice  consisting  of,  alternately,  Alo.25Gao.75 As  and  GaAs  layers.  For  the  purpose  of 
this  paper,  we  have  chosen  the  thickness  of  the  Alo.25Gao.75As  layer  as  6  nm  and  because 
of  our  familiarity  with  it,  the  Matlab  programming  language  has  been  used. 

A.  BANDSTRUCTURE  FUNCTION 

Before  presenting  the  objective  of  this  paper,  i.e.,  plotting  the  Wannier  functions 
of  the  GaAs/Al0.25Gao.75As  superlattice,  we  wanted  to  make  sure  that  the  eigenvalue 
matrix  equation  (Eqn.  32)  that  we  formed  as  the  discretization  of  the  superlattice 
Schrodinger  equation  (Eqn.  27)  is  accurate.  The  easiest  and  the  most  effective  way  to  do 
so  would  be  to  compare  the  eigenvalues  of  the  matrix  equation  (Band  structure  function) 
to  another  method  (which  has  already  proven  itself  to  be  correct)  of  finding  the  band 
structure  function  of  a  superlattice  structure. 

We  first  created  a  subroutine  in  Appendix  A  for  the  evaluation  of  the  eigenvalues 
e n  (k)  of  the  matrix  equation  (Eqn.  32)  which  gave  us  the  bandstructure  function  of  the 
superlattice.  Next  in  the  second  subroutine,  we  used  G.  Bastard’s  analytic  superlattice 
dispersion  relation  equation  [Ref.  24] 

\(  k  k  \ 

Cos(kd)  =  Cos(kGdG )Cos{kAdA) - ~+~j~  Sin{kcdc )SinikAdA )  (33) 

kA  J 


17 


where 


k  -k  = 

kG  kGoAs 


2  Sn(k)mGaAs 


k  -k 

kA  kAIGoAs 


2en(k)mAlGaAs 


=  ,/AZs„  + 


ti2 


^G  ~  d GaAs 


>  di  ~  d 


A  ~  u  At, GaAs 


(34) 


to  compare  with  the  eigenvalues  of  the  matrix  equation.  Equation  33  and  34  are  implicit 
equations  relating  the  energy  values  en  (k)  to  the  wave  vector  k. 


Figure  2  shows  the  bandstructure  function  obtained  by  evaluating  the  eigenvalues 
of  the  matrix  equation  (Eqn.  32),  while  Figure  3  is  the  solution  of  the  Equation  34 
satisfying  the  dispersion  relation  (Eqn.  33).  By  inspection,  it  is  clearly  seen  that  the  two 
bandstructure  figures  are  in  good  agreement  with  each  other.  That  is  to  say,  the 
eigenvalue  matrix  equation  formed  as  the  discretization  of  the  superlattice  effective-mass 
Schrodinger  equation  (Eqn.  27)  is  accurate. 


BANDSTRUCTURE  OF  GaAs/AIGaAs(AI%25)  SUPERLATTICE 


Figure  2:  Bandstructure  function  £n  (k)  obtained  by  evaluating  the 
eigenvalues  of  the  matrix  equation  (Eqn.  32)  where  A Ec  =  E*!GaAs  —  EGaAs . 
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V  Oo 


BANDSTRUCTURE  OF  GaAs/AIGaAs  SUPERLATT1CE  (Eqn.  33) 


Figure  3:  Bandstructure  function  £n(k)  obtained  by  solving  the  G.  Bastard’s 
superlattice  dispersion  relation  equation  (Eqn.  33  and  34)  where 

_  j?  AlGaAs  _ jjGaAs 


B.  WANNIER  FUNCTIONS 

After  confirming  the  accuracy  of  the  numerical  modeling  of  the  superlattice 
Schrodinger  equation  (by  comparing  with  the  G.  Bastard’s  superlattice  dispersion 
relation),  the  computer  code  in  Appendix  B  implements  an  algorithm  which  calculates 
and  plots  the  Wannier  functions  of  the  GaAs/Alo.25Gao.75As  superlattice  structure.  At 
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this  phase  of  our  numerical  modelling,  we  deliberately  skip  the  normalization  of  the 
Bloch  waves.  There  are  three  reasons  for  this:  First  of  all,  the  Bloch  functions  can  easily 
be  calculated  and  have  been  obtained  by  solving  the  Schrodinger  equation  for  each  value 
of  k  and  n  [Ref.  25].  Secondly,  we  are  not  concerned  with  the  normalization  of  the  Bloch 
functions,  since  having  arbitrary  amplitudes  for  the  Wannier  functions  of  any  band  index 
(n)  is  immaterial  for  the  purpose  of  this  paper,  i.e.,  plotting  and  pointing  out  the 
properties  of  the  Wannier  functions.  The  third  reason  which  is  the  limitation  of  our 
modelling  is  that  normalizing  each  Bloch  function  of  every  superlattice  unit  cell  would  be 
formidable  to  carry  out  numerically  if  you  discretized  300  mesh  points  for  each 
superlattice  site  (zj)  and  100  mesh  points  in  the  first  BZ  (k)  which  brings  out  30,000 
Bloch  functions  just  for  the  first  band  index  (n=  1)  to  deal  with.  We  are  going  to  point  out 
this  computational  difficulty  in  the  conclusion  and  final  remarks  subchapter. 

In  Chapter  II,  we  have  with  inclusive  detail  proved  that  for  each  band  index  («) 
there  exists  only  one  Wannier  function.  Since  for  most  of  the  optical  and  electronic 
applications  the  first  two  bands  of  the  bandstructure  function  are  utilized,  in  the  computer 
code  (Appendix  B)  we  numerically  calculated  the  Wannier  functions  of  those  bands,  i.e., 
n=  1  and  n= 2.  Favorably,  plotting  the  Wannier  functions  for  these  two  bands  suffices  for 
the  purpose  of  this  paper.  Furthermore,  we  developed  the  computer  code  in  such  a  way 
that  by  changing  the  band  index  entry  in  the  program  we  can  obtain  the  Wannier 
functions  of  any  band  index  that  we  are  interested  in. 

Figure  4  shows  the  Wannier  function  of  the  first  band  (aj(z))  extending  to  the 
right  neighboring  superlattice  site.  By  inspection,  it  is  clearly  recognized  that  the 
Wannier  function  (aj(z))  is  localized  and  centered  at  the  superlattice  site  which  is  shown 
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WANNIER  FUNCTION  OF  GaAs/AIGaAs(AI%25)  SUPERLATTICE 


Figure  4  :  The  Wannier  function  of  GaAs/Al^  Gao.75  As  superlattice  for  the  first 
band  index  (arfz)). 


with  a  thick  line  in  Figure  4.  The  asymptotic  behavior  is  also  observed  easily  as  the 
W annier  function  decays  down  in  the  neighboring  superlattice  unit  cell  since  as 

increases  (Eqn.  9),  the  Bloch  functions  destroy  each  other  by  interference  which  we  had 
proved  thoroughly  in  Chapter  II  (  See  the  discussion  after  Equation  9). 

We  shall  now  show  that  the  Wannier  functions  are  real.  That  is  to  say,  if  we  tried 
to  plot  the  imaginary  part  of  the  numerically  calculated  Wannier  function  of  the  first  band 
(or  any  band),  we  should  get  an  empty  figure.  However,  since  the  numerical  approach 
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(finite  difference  method)  that  we  used  to  solve  the  superlattice  Schrodinger  equation 
(Eqn.  27)  is  just  an  approximation,  we  would  expect  some  amplitude  as  error.  Figure  5 


0  0.5  1  1.5  2  2.5  3 


z ,  (nm) 

Figure  5:  Imaginary  part  of  the  numerically  calculated  Wannier  function  of 
GaAs/Alo.25  Gao.75  As  superlattice  for  the  first  band  index  (n=l). 

shows  the  imaginary  part  of  the  Wannier  function  of  the  first  band  (arfz)).  Indeed,  if  we 
note  the  order  of  amplitude  which  is  around  10‘15,  we  can  confidently  state  that  the 
Wannier  functions  for  theGaAs/Alo.25Gao.75As  superlattice  is  real. 
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Figure  6  depicts  the  Wannier  function  of  the  superlattice  (SL)  structure  for  the 
band  index  n-2  (a2(z))-  The  localised  nature  and  the  asymptotic  behavior  are  obviously 
recognized  as  we  did  for  the  ai(z).  Again,  the  Wannier  function  a2(z)  is  centered  at  the 


WANNIER  FUNCTION  OF  GaAs/AIGaAs(AI%25)  SUPERLATTICE 


Figure  6:  The  Wannier  function  of  GaAs/Al0.25  Gao.75  As  superlattice  for  the 
second  band  index  (a2(z)). 

same  superlattice  unit  cell  as  we  have  verified  for  at(z)  in  Figure  4.  The  important  aspect 
of  Figure  6  is  the  anti-symmetry  of  the  Wannier  function  for  the  band  index  2  ( a2(z. )) 
about  z=0  where  we  can  easily  see  the  symmetry  of  the  Wannier  function  for  the  band 
index  n=l(at(z))  in  Figure  4.  W.  Kohn  has  given  a  comprehensive  discussion  about  the 
symmetry  properties  of  the  Bloch  and  the  Wannier  functions  [Ref.  25],  Basically  the 
(anti-)  symmetry  property  of  the  Wannier  functions  about  z= 0  can  be  accounted  for  by 
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Amplitude  (unitless) 


the  parity  of  the  Bloch  functions  which  could  be  easily  perceived  if  we  revisited  the  basic 
quantum  mechanics  problem  of  square  well. 

The  reality  of  the  Wannier  functions  can  again  be  shown  if  we  plot  the  imaginary 
part  of  a2(z)  which  is  resulted  from  our  numerical  approach  (Fig.  7).  Note  the 


10-15  WANNIER  FUNCTION  (Imaginary)  OF  GaAs/AIGaAs(AI%25)  Si- 


Figure  7:  Imaginary  part  of  the  numerically  calculated  Wannier  function  of 
GaAs/Alo.25  Gao.75  As  superlattice  for  the  first  band  index  (n= 2). 


amplitude  of  the  numerically  calculated  imaginary  part  of  a2(z)  which  is  around  10'15. 

C.  CONCLUSION  AND  FINAL  REMARKS 

In  this  thesis,  we  have  developed  an  algorithm  to  numerically  calculate  the 
Wannier  functions  of  a  GaAs/Alo.25  Gao.75  As  superlattice  structure.  The  algorithm  has 
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been  built  in  two  seperate  steps :  First  we  wanted  to  show  the  credibility  of  our 
numerical  modelling  of  the  superlattice  effective-mass  Schrodinger  equation  (Eqn.  27). 
We  did  that  by  comparing  the  bandstructure  function  modelled  by  our  method  with  the 
one  obtained  by  using  the  G.  Bastard’s  superlattice  dispersion  relation  (Appendix  A). 
We  proved  that  our  modelling  with  the  finite  difference  form  is  a  very  good 
approximation  for  the  discretization  of  the  superlattice  effective-mass  Schrodinger 
equation  (Figures  2  and  3).  In  the  next  step,  after  confirming  that  our  numerical 
modelling  is  reliable,  we  created  a  routine  (Appendix  B)  for  solving  the  eigenvalue 
matrix  equation  (Eqn.  32)  for  the  Wannier  functions  of  the  GaAs/Alo.25  Gao^sAs. 

The  results  we  obtained  are  in  good  agreement  with  the  theory  of  the  Wannier 
functions  which  we  have  discussed  in  detail  in  Chapter  II.  Although,  in  the  past  there 
were  no  attempts  to  calculate  these  functions  for  bulk  crystals,  we  have  demonstrated 
that  the  Wannier  functions  can  be  calculated  even  for  a  complex  superlattice  structure. 
By  plotting  the  numerical  simulation  of  the  Wannier  functions  for  the  first  couple  of 
allowed  bands  ( n=l  and  n= 2)  of  the  superlattice  structure  ,  we  displayed  the  peculiar 
properties  of  these  functions: 

a.  Wannier  functions  are  real  (Figures  5  and  7). 

b.  They  fall  off  exponentially  (Figures  4  and  6). 

c.  They  are  either  symmetric  or  anti-symmetric  about  z=0  (Figures  4  and  6). 

Even  though  the  results  we  reach  conform  with  the  theory,  the  numerical 

approach  we  have  utilized  create  two  shortcomings  that  we  have  to  deal  with.  The  first 
one  which  we  pointed  out  before  is  the  normalization  problem  of  the  Bloch  functions. 
Since  for  every  band  index  n  and  for  each  value  of  k  there  is  a  Bloch  state  ( y/nk  (z) ), 
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normalization  for  a  continuous  k  values  is  a  formidable  task  to  do,  e.g.,  just  for  the  first 
band,  the  computer  code  in  Appendix  B  has  simulated  16800  Bloch  functions  for  discrete 
values  of  k.  Besides,  we  have  to  note  that  we  just  extend  the  Wannier  functions  to  the 
next  two  superlattice  sites  to  show  the  asymptotic  behavior.  But,  as  we  have  stated 
before,  for  the  purpose  of  this  paper  it  is  unnecessary  to  normalize  the  Bloch  states.  The 
second  problem  we  encounter  is  about  the  effective  masses.  In  order  to  show  all  the 
aspects  of  the  Wannier  functions  we  have  had  to  use  a  uniform  effective  mass  throughout 
the  simulation  distance.  The  reason  for  that  is  the  same  as  the  first  one:  We  could  not 
discretize  the  superlattice  site  small  enough  to  show  the  kink  at  each  interface  the 
Wannier  functions  should  have  had  that  would  be  due  to  the  different  effective  masses  at 
the  interface.  We  used  the  smallest  mesh  distance  that  we  could,  but  all  we  obtained  was 
discontinuities  at  the  interfaces. 

Finally,  we  want  to  call  attention  to  some  of  the  future  studies  that  we  think  will 
support  the  theory  of  the  Wannier  functions: 

a.  We  have  used  Bloch  states  in  the  eigenvalue  matrix  equation  to  obtain  the 
bandstructure  function  of  the  GaAs/Alo.25  GaojsAs  superlattice  structure.  An 
algorithm  can  be  developed  to  use  the  Wannier  functions  in  the  calculation  of 
the  energy  bandstructure  functions  of  the  superlattices. 

b.  AlxGai.xAs  bulk  alloy  is  a  direct-bandgap  semiconductor  if  x<0.45  [Ref.  26] 
There  is  a  crossover  from  direct  to  indirect-bandgap  structure  at  about  x=0.45. 
The  Wannier  functions  of  such  a  structure  can  be  studied  in  the  future. 
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c.  The  GaAs/Alo.25  Gao  75AS  superlattice  structure  is  a  type  1  superlattice.  The 
Wannier  functions  of  type  2  and  type  3  superlattices  can  be  studied  in  the 
future. 

d.  If  Wannier  functions  of  all  types  of  superlattice  structures  are  obtained,  maybe 
an  electronic  transport  theory  can  be  developed  based  on  the  Wannier 
functions. 
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APPENDIX  A:  COMPUTER  CODE  IN  MATLAB  PROGRAMMING 
LANGUAGE  FOR  NUMERICAL  CALCULATION  OF  THE 
BANDSTRUCTURE  FUNCTION 


format  long  g 


%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 
%%%%%%%%%%%%%%%%%%%%%%%  INPUT  PARAMETERS  %%%%%%%%%%%%%%%%%%%%%%%%% 
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 

EgG= input ( 'Energy  gap  for  GaAs  (eV) ' ) ; 
hh=input ( 'barier  thickness  (nm)  ' ) ; 

LL= input ( 'well  thickness  (nm) ' ) ; 

x= input ( '  A1  concentration  (x)  ' )  ; 


%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 
%%%%%%%%%%%%%%%  DEFINING  THE  CONSTANTS  %%%%%%%%%%%%%%%%% 
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 

hbar=6 . 5821e-16; 
ml=0 . 067 ; 

me=(0.511e6)/9el6; 
m2  =  (ml+0 . 0838*x) ; 

EgA=EgG+1.247*x; 

Ec= (0.65* (EgA-EgG) ) ;  %CONDUCTION  BAND  OFFSET 

h=hh/ 100 00 00000 ; 

L=LL/1000000000 ; 
d=h+L ; 

D= (h/2 ) +L; 

s=3 . 14/d; 

in=100; 

dz=d/in; 

dq=f loor (s/50)  ; 

w=-  (hbar*2)  /  (2*me*  (dz^2)  )  ; 


%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 
%%%%%%%%%%%%%%%%%%%%%  HAMILTONIAN  (H  )  %%%%%%%%%%%%%%%%%% 
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 

for  1=1: in 

z (1) =dz*l; 

if  z (1) >h/2  &  z (1) <D 
V(l) =0 ; 

HH  (1)  =V (1)  -  (2 /ml)  *w; 
m  ( 1 ) =ml ; 

Ww  ( 1 )  =w/ml ; 

else 

V  ( 1 ) =Ec ; 

HH(l)=V(l)-(2/m2)*w; 

Ww(l)  =w/m2 ; 

end 

end 

Ww( floor (in/ 6) ) =  (2*w) / (ml+m2) ;Ww( floor (in*5/6) ) = (2*w) / (ml+m2) ; 

HHl=diag (HH) ? 

W=Ww ( 1 : ( length (Ww) -1 ) ) ; 

Wl=diag (W,l) ; 
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W2=diag(W,-l) ; 
HH2= (HH1+W1+W2 ) ; 


%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 
%%%%%%%%%%%%%%%%  FIRST  BROULINE  ZONE  (BZ)  %%%%%%%%%%%% 

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 


for  k= 1:50; 


end 


ql (k) = (k) *dq; 

Ql= (w/m2) * (cos (ql (k) *d) -i*sin(ql (k) *d) ) ; 
Q2= (w/m2) * (cos (ql (k) *d)+i*sin(ql (k) *d) ) ; 
Qll=diag(Ql, (in-1) ) ;Q22=diag (Q2 , - (in-1) ) ; 


H= (HH2+Q11+Q22 ) ; 

kk=eig(H) ; 
ee=sort(kk) ; 


%  HAMILTONIAN  WITH  BOUNDARY  CONDITIONS 
%  ENERGY  EIGENVALUES 


EE1 (k) =ee (1) 
EE2 (k) =ee (2 ) 
EE3 (k) =ee (3 ) 


%  BAND  INDEX  n=l 
%  BAND  INDEX  n=2 
%  BAND  INDEX  n=3 


qql=-l*f liplr (ql) ; 
q= [qqi  ql] ; 

EEll=fliplr (EE1) ;E1= [EE11  EEl]  ; 
EE22=f liplr (EE2 ) ;  E2=[EE22  EE2 ] ; 


EE33=fliplr (EE3) ;E3=[EE33  EE3] ; 

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 
%%%%%%%  PLOTTING  BANDSTRUCTURE  FUNCTION  %%%%%%%%%%%%%% 
%%%%%%%  WITH  EIGENVALUE  MATRIX  EQUATION  %%%%%%%%%%%%%% 
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 

figure (1) (plot(q,El, ' k' , q, E2 , ' k' , q, E3 ,  'k' )  ; 
axis ( [-s  s  0  . 45] ) 

title ( 'BANDSTRUCTURE  OF  GaAs/AlGaAs (Al%25)  SUPERLATTICE ' ) 
xlabel('k  (first  BZ) ,  1/m') 
ylabel ( 'Energy,  eV' ) 


%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 


%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 
%%%%%%%%  G.  BASTARD'S  DISPERSION  RELATION  %%%%%%%%%%% 
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 


mml=ml* ( 0 . 511e6) /9el6 ; 

mm2= (ml+0 . 0838*x) * (0 . 511e6) /9el6 ; 

wl=sqrt (2*mml) /hbar; 
w2=sqrt (2*mm2) /hbar; 

dq=1000000 ; 
dE=0. 000002; 

d=h+L;dz=d/ 600;D=h/2+L; 
it=ceil (s/dq) ; 

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 
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%%%%%%%%%%%%%%%%%  DISPERSION  RELATION  %%%%%%%%%%%% 
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 


n=l : 1 : 200000; 

E= . 05+n . *dE; 

pl=(cos( (sqrt(E) ) ,*wl*L) ) . * (cosh ( (sqrt (Ec-E)  )  ,*w2*h) ) ; 
p2=l/2* ( (-sqrt ( (Ec-E) ./E) . *sqrt (m2/ml) )  +  (1 . / (sqrt ( (Ec- 
E) . /E) . *sqrt (m2/ml) ) ) ) ; 

p3= (sin ( (sqrt (E) ) .*wl*L) ) .  * (sinh( (sqrt (Ec-E) ) . *w2*h)  )  ; 
p=pl-p2 . *p3 ; 

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 
%%%%%%%%%%%%%%%%  FIRST  BROULINE  ZONE  (BZ)  %%%%%%%%%%%% 
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 


for  v=l : it 
q=v*dq; 

g(v)=cos(q*d) ; 

[yy  hh]=min(abs(g(v)-p) )  ; 

Ee (v) = . 05+hh*dE; 
qq(v)=q; 

end 

f f= (-1) *fliplr (qq) ; 

k=[ff  qq]  %  FIRST  BRUOILLOIN  ZONE 

ee=f liplr (Ee) ; 

EE= [ee  Ee] ;  %  E(k) 

figure(2) ,plot(k,EE, 'b. ' ) 
axis([-s  s  0  0.45]) 

title ([' BANDSTRUCTURE  OF  GaAs/AlGaAs (x= . 25)  SUPERLATTICE  (Bastard)']) 
xlabel('k  (first  BZ) ,  1/m' ) 
ylabel ( 'Energy  (eV) ' ) 
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APPENDIX  B:  COMPUTER  CODE  IN  MATLAB  PROGRAMMING 
LANGUAGE  FOR  NUMERICAL  CALCULATION  OF  THE 
WANNIER  FUNCTIONS 


format  long  g 

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 
%%%%%%%%%%%%%%%%%%%%%%%  INPUT  PARAMETERS  %%%%%%%%%%%%%%%%%%%%%%%%% 
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 

EgG=input ( 'Energy  gap  for  GaAs  (eV)  ' )  ; 
hh= input ( 'barier  thickness  (nm) ' ) ; 

LL= input ( 'well  thickness  (nm) ' ) ; 
x=input ( ' Al  concentration  (x) ' ) ; 

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 
%%%%%%%%%%%%%%%%%%  DEFINING  THE  CONSTATNT S  %%%%%%%%%%%%%%%%%%%%%%%% 
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 


hbar=6 . 5821e-16 ; 
ml=0 . 067 ; 

me= ( 0 . 511e6 ) /9el6; 
m2=  (ml+0 . 0836*x) ; 

EgA=EgG+l . 247*x;  %AlGaAs  bandgap  for  given (x) 

Ec= (0.65* (EgA-EgG) ) ;  %Conduction  band  offset 

h=hh/ 1000000000  ? 

L=LL/1000000000; 
d= (h) +L; 

D= (h/2 ) +L;  DH=D+h ;  DHL=DH+L  ; DHLH=DHL+ (h) ;  DHLHL=DHLH+L ; 

s=3 . 14/d; 

ind=210? 

dz=d/70 ; 

inq=40 ; 

dq=floor (s/inq) ; 

w=-  (hbar'N2)  /  (2*me*  (dz/v2 ) )  ; 


%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%.%%%%%%%%%%%%%%%%%%%% 
%%%%%%%%%%%%%%%%%%%%%  HAMILTONIAN  (H)  %%%%%%%%%%%%%%%%%%%%%%%% 

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 


for  l=l:ind 
z (1) =dz*l; 

if  z (1) > (h/2 )  &  z ( 1 ) <D 
V (1) =0 ; 

HH  (1)  =V(1)  -  (2/ml)  *w; 
m ( 1 ) =ml ; 

Ww ( 1 ) =w/ml ; 

elseif  z ( 1 ) >D  &z(l)<DH 
V (1) =Ec ; 

HH  ( 1 )  =V  ( 1 )  -  ( 2  / ml )  * w  ; 

Ww  ( 1 )  =w/ml ; 

elseif  z ( 1 ) >DH  &  z(l)<DHL 
V  ( 1 )  =  0  ; 

HH  ( 1 )  =V  ( 1 )  -  ( 2  /ml )  *w; 
m  ( 1 ) =ml ; 

Ww  ( 1 )  =w/ml ; 
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elseif  z(l)>DHL  &  z(l)<DHLH 
V ( 1 ) =Ec ? 

HH  ( 1 ) =V ( 1 ) -  ( 2  /ml )  *w; 

Ww ( 1 ) =w/ml ; 

elseif  z ( 1 ) >DHLH  &  z(l)<DHLHL 
V  ( 1 )  =  0 ; 

HH ( 1 )  =V ( 1 )  -  {2 /ml)  *w; 
m  ( 1 ) =ml  ; 

Ww ( 1 ) =w/ml ; 

else  % z { 1 ) >DLHL  %&  z{l)<DLHLD 
V (1) =Ec ; 

HH  ( 1 )  =V  ( 1 )  -  ( 2  /  ml )  *w ; 

Ww ( 1 ) =w/ml ; 

end 

end 


%Ww  { floor  ( ind/18 ) )  =  { 2  *w)  /  (ml+m2 )  ; 
%Ww(f loor (ind*7/18) )  =  (2*w) / (ml+m2) ; 
%Ww{  floor  (ind*  11/18) )  -  (2*w)  /  (ml+m2 )  ; 
%Ww  ( floor  ( ind*  13 / 18 )  )  =  ( 2 *w)  /  (ml+m2 )  ? 
%Ww  (floor  ( ind*  17/18)  )  =  (2*w)  /  (ml+m2)  ; 


HHl=diag (HH) ; 

W=Ww ( 1 : (length (Ww) -1) )  ; 

Wl=diag (W, 1) ? 

W2=diag(W,-l); 

HH2= (HH1+W1+W2 ) ? 

nl=zeros (ind, 1) ;n2=zeros (ind, 1) ;n3=zeros (ind,  1) ; 

for  k=l:2*inq; 
if  k>inq 

ql (k-inq) = (k-inq) *dq; 

Ql= (w/ml) * (cos (ql (k-inq) *d) -i*sin (ql (k-inq) *d) ) ; 
Q2= (w/ml) * (cos (ql (k-inq) *d) +i*sin (ql (k-inq) *d) ) ; 

else 

ql (k) =-k*dq; 

Ql= (w/ml) * (cos (ql (k) *d) -i*sin(ql (k) *d) ) ; 

Q2= (w/ml) * (cos (ql (k) *d) +i*sin(ql (k) *d) ) ; 
end 

Qll=diag(Ql/ (ind-1) ) ;Q22=diag (Q2 , - (ind-1 ) ) ? 

H= (HH2+Q11+Q22) ; 

[jj  kk] =eig (H) ; 
ee=sort (diag (kk) ) ; 

[fl  gl ] =f ind ( ee ( 1 ) ==kk) ; 

[ f 2  g2] =f ind (ee (2) ==kk) ; 

E f 3  g3 3 =f ind (ee (3 ) ==kk) ; 

nl=nl+ ( j j ( : , fl) ) ; 
n2=n2+ ( j j ( : , f2) ) ; 

end 


%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 
%%%%%%%%%%%%%%  PLOTTING  OF  THE  WANNIER  FUNCTIONS  %%%%%%%%%%%%%%%%%%% 
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%^%^'^^^^%%%%%%%%%%%% 


g=l : ind; 

figured)  rplot  {  (g.*dz)  .*100000000^1,  'k' ) 

title ( 'WANNIER  FUNCTION  OF  GaAs/AlGaAs (Al%25)  SUPERLATTICE ' ) 
xlabel('z  ,  nm  ' ) 
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ylabel ( 'Amplitude  (arbitrary) ' ) 

figure (2 ) ,plot( (g.*dz) . *100000000, n2, 'k' ) 

title ('WANNIER  FUNCTION  OF  GaAs/AlGaAs (Al%25)  SUPERLATTXCE  ') 

xlabel  ( '  z  ,  run  ' ) 

ylabel ( 'Amplitude  (arbitrary) ' ) 
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